Particle methods are known as simulation methods for representing a motion of a continuum such as a fluid, an elastic body or the like, and a powder and granular assembly by using a particle swarm. The particle methods include Moving Particle Semi-implicit (MPS) method, Smoothed Particle Hydrodynamics (SPH) method, and the like. These particle methods can more easily handle problems including a moving boundary of a free-surface flow or the like in comparison with simulation methods such as a finite difference method, a finite element method, and the like, which use a numerical mesh. Therefore, the particle methods have been increasingly used in various fields in recent years.
The particle methods are methods developed to solve a problem including a moving boundary of a free-surface flow, or the like, and a term that is called an artificial viscosity term and has an effect of suppressing a motion between particles is added to a calculation formula in order to stably calculate a motion of a particle.
For example, with the SPH method, a momentum conservation law and amass conservation law are discretized as follows.
                                          ⅆ                          v              a                                            ⅆ            t                          =                  -                                    ∑              b                        ⁢                                                  ⁢                                                            m                  b                                [                                                      (                                                                                            P                          b                                                +                                                  P                          a                                                                                                                      ρ                          b                                                ⁢                                                  ρ                          a                                                                                      )                                    +                                      Π                    ab                                                  ⁢                                                                  ]                            ⁢                                                ∂                                      W                    ⁡                                          (                                                                                                                                                            r                              a                                                        -                                                          r                              b                                                                                                                                ,                        h                                            )                                                                                        ∂                                      r                    a                                                                                                          (        1        )                                                      ⅆ                          ρ              a                                            ⅆ            t                          =                              ∑            b                    ⁢                                          ⁢                                    m              b                        ⁢                                          ρ                a                                            ρ                b                                      ⁢                                          (                                                      v                    a                                    -                                      v                    b                                                  )                            ·                                                ∂                                      W                    ⁡                                          (                                                                                                                                                            r                              a                                                        -                                                          r                              b                                                                                                                                ,                        h                                            )                                                                                        ∂                                      r                    a                                                                                                          (        2        )            
In the formulas (1) and (2), a subscript a indicates an a-th particle (particle a), and a subscript b indicates a b-th particle (particle b). ra, va, ρa, Pa, and ma respectively indicate a position vector, a velocity vector, a density, a pressure, and a mass of the particle a. Similar notations are used also for the particle b.
  ∑  bon the right side of the formulas (1) and (2) indicates a total sum with respect to all particles other than the particle a.
W is a kernel function, and used to configure a continuous field from a distribution of particles. The following cubic spline function or the like is used as the kernel function W in many cases.
                              W          ⁡                      (                          r              ,              h                        )                          =                  {                                                                                          (                                          1                      -                                              1.5                        ⁢                                                                              (                                                          r                              h                                                        )                                                    2                                                                    +                                              0.75                        ⁢                                                                              (                                                          r                              h                                                        )                                                    3                                                                                      )                                    /                  ϕ                                                                              0                  ≤                                      r                    h                                    <                  1                                                                                                      0.25                  ⁢                                                                                    (                                                  2                          -                                                      r                            h                                                                          )                                            3                                        /                    ϕ                                                                                                1                  ≤                                      r                    h                                    <                  2                                                                                    0                                                              2                  ≤                                      r                    h                                                                                                          (        3        )            
h indicates a radius of an influence between particles, and a value as large as twice to three times of a mean interval between particles in an initial state is used in many cases. φ is a normalization factor adjusted so that an integrated value of the entire space of the kernel function W results in 1. This factor is decided to be 0.7 πh2 in a case of a two dimensional space, and to be πh3 in a case of a three dimensional space.
Πab in the second term on the right side of the formula (1) is an artificial viscosity term that works between the particles a and b. For example, the following form is adopted.
                              Π          ab                =                  {                                                                                                                                        -                        α                                            ⁢                                                                                          ⁢                      c                      ⁢                                                                                          ⁢                                              μ                        ab                                                              +                                          βμ                      ab                      2                                                                            ρ                    ab                                                                                                                                          (                                                                        v                          a                                                -                                                  v                          b                                                                    )                                        ·                                          (                                                                        r                          a                                                -                                                  r                          b                                                                    )                                                        <                  0                                                                                    0                                                                                                        (                                                                        v                          a                                                -                                                  v                          b                                                                    )                                        ·                                          (                                                                        r                          a                                                -                                                  r                          b                                                                    )                                                        ≥                  0                                                                                        (        4        )                                          μ          ab                =                                            h              ⁡                              (                                                      r                    a                                    -                                      r                    b                                                  )                                      ·                          (                                                v                  a                                -                                  v                  b                                            )                                                                                                                                r                    a                                    -                                      r                    b                                                                              2                        +                          η              2                                                          (        5        )                                          ρ          ab                =                                            ρ              a                        +                          ρ              b                                2                                    (        6        )            
α, β, and η are constants, and c indicates a sound velocity.
The formula (1) is an equation that represents a momentum conservation law of a fluid. The first term on the right side of the formula (1) represents a force applied by a pressure of the fluid, and has an effect such that particles press against each other. The artificial viscosity term Πab has an effect of damping a momentum as a brake applied when particles get close each other, and of preventing the particles from moving through. The formula (2) is an equation that represents a mass conservation law of the fluid, and has an effect such that a density increases when particles get close each other or decreases when particles draw apart.
However, the forms of the formulas (1) to (3) using the artificial viscosity term make the effect that a motion is damped by the artificial viscosity term significant. Accordingly, also Godunov SPH method and the like that are intended to suppress damping of a motion by the artificial viscosity term to a minimum and use a solution to a Riemann problem have been devised.
Additionally, a particle-type fluid simulation method for reproducing wetness and repelling between a fluid and a surface of a solid, and a data processing method for efficiently executing a computational process for data of particles present in meshed space have been devised.
Patent Document 1: Japanese Laid-open Patent Publication No. 10-185755
Patent Document 2: International Publication Pamphlet No. WO 2008/020634
Non-patent Document 1: S. Koshizuka, H. Tamako and Y. Oka, “A Particle Method for Incompressible Viscous Flow with Fluid Fragmentation”, Computational Fluid Dynamics Journal, vol. 4, no. 1, p. 29-46, April 1995
Non-patent Document 2: J. J. Monaghan, “Smoothed Particle Hydrodynamics”, Annual Review of Astronomy and Astrophysics, vol. 30, p. 543-574, 1992
Non-patent Document 3: S. Inutsuka, “Reformulation of Smoothed Particle Hydrodynamics with Riemann Solver”, Journal of Computational Physics 179, p. 238-267, 2002